Tailoring Magnetic Doping in the Topological Insulator Bi 2 Se 3 
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We theoretically investigate the possibility of establishing ferromagnetism in the topological insulator Bi2Se3 
via magnetic doping of 3d transition metal elements. The formation energies, charge states, band structures, 
and magnetic properties of doped Bi2Se3 are studied using first-principles calculations within density functional 
theory. Our results show that Bi substitutional sites are energetically more favorable than interstitial sites for 
single impurities. Detailed electronic structure analysis reveals that Cr and Fe doped materials are still insulating 
in the bulk but the intrinsic band gap of Bi 2 Se3 is substantially reduced due to the strong hybridization between 
the d states of the dopants and the p states of the neighboring Se atoms. The calculated magnetic coupling 
suggests that Cr doped Bi2Se3 is possible to be both ferromagnetic and insulating, while Fe doped Bi2Se3 tends 
to be weakly antiferromagnetic. 
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The recent discovery of topological insulators (TIs) in nar- 
row band gap semiconductors has generated great interest 
in both condensed matter physics and materials science [1- 
7]. The strong spin-orbit coupling (SOC) in these materials, 
which is responsible for the nontrivial band topology, also 
highlights their application potential in semiconductor spin- 
tronics [8]. In particular, the intricate interplay between topo- 
logical order and ferromagnetism is expected to give rise to a 
variety of unconventional spintronic effects that may lead to 
entirely new device paradigms [9-12]. For example, a weak 
magnetic perturbation can open up an energy gap in the sur- 
face spectrum of a TI, resulting in the emergence of massive 
Dirac fermions [13]. In a TI thin film, long range magnetic 
order can also be established, which potentially could realize 
the long-sought quantum anomalous Hall state [14-16]. 

A natural strategy to introduce magnetism in TIs is via 
magnetic doping, similar to diluted magnetic semiconduc- 
tors [17]. The current research mainly focuses on the Bi 2 Se3 
family [13-15, 18-20]. However, existing experimental evi- 
dences have shown that the solubilities of magnetic impurities 
are normally very limited in Bi2Se3 family compounds [13, 
18, 21, 22], magnetic impurities tend to form clusters instead 
of being individually distributed in the materials, and dopants 
have uncertain valency [23]. On the magnetic ordering, ferro- 
magnetism was reported in Mn and Fe doped Bi2Te3 [22, 24- 
27] and V, Cr and Mn doped Sb 2 Te 3 [21, 22, 25], while a 
spin glass state was found in Mn doped Bi2Se3 [22]. For Cr 
doped Bi 2 Se3, both antiferromagnetism [28] and ferromag- 
netism [29] were observed, and the observations in Fe doped 
Bi 2 Se3 were also rather controversial [24, 28, 30, 31]. From 
a growth point of view, it is imperative to choose proper mag- 
netic elements and be able to precisely control the magnetic 
impurity distribution inside the host materials [13]. This re- 
lies on comprehensive knowledge of the kinetic and energetic 



behaviors of magnetic impurities in TIs and consequently the 
electronic and magnetic properties. Yet, such knowledge is 
still lacking in the literature. 

In this Letter, we investigate the feasibility of tailing Bi 2 Se3 
to a ferromagnetic insulator via doping 3d transition metal 
(TM) elements, including V, Cr, Mn, and Fe, using first- 
principles calculations within density functional theory. We 
find that Bi substitutional sites are energetically more favor- 
able than interstitial sites for single impurities of these ele- 
ments, and the optimal growth conditions for TM doping are 
presented. The calculated band structures reveal that only Cr 
and Fe doped materials are insulating in the bulk but the in- 
trinsic band gap of Bi 2 Se3 is substantially reduced due to the 
strong hybridization between the d states of the dopants and 
the p states of the neighboring Se atoms. Further investiga- 
tion on the magnetic coupling between the dopants clarifies 
some experimental discrepancies and suggests that Cr doped 
Bi 2 Se3 is likely to be ferromagnetic, while Fe doped material 
tends to be weakly antiferromagnetic. 

Our first-principles density functional theory calcula- 
tions are performed using the pseudopotential plane-wave 
method with projected augmented wave [32] potentials and 
Perdew-Burke-Ernzerhof-type generalized gradient approxi- 
mation (GGA) [33] for exchange-correlation functional, as 
implemented in the Vienna ab initio simulation package 
(VASP) [34, 35]. A plane- wave energy cutoff of 300 eV is 
consistently used in all the calculations. Unless mentioned 
otherwise, SOC [36] is taken into account in all calcula- 
tions. Crystalline Bi 2 Se3 has a rhombohedral structure and 
its unit cell is composed of three weakly coupled quintuple 
layers (QLs). To investigate the behaviors of isolated mag- 
netic dopants in bulk Bi 2 Se3, we employ a 2 x 2 x 1 supercell 
containing 24 Bi and 36 Se atoms, as shown in Fig. 1 (a). A 
gamma-centered 7x7x2 mesh of special k points is adopted 
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FIG. 1: (color online), (a) Illustration of a 2 x 2 x 1 supercell for 
modeling a single dopant in bulk Bi2Se3. (b) Calculated formation 
energies of the most stable configurations of single V, Cr, Mn, and 
Fe impurities doped Bi 2 Se3 as a function of the host element chem- 
ical potentials. For comparison, the values of dominant intrinsic de- 
fects, selenium vacancies Vsei, antisite defects Bisei and Sesi, are 
shown in dashed lines. 



for integrations over the Brillouin zone. In calculating the for- 
mation energies [37, 38] of charged defects, a uniform back- 
ground charge is introduced to keep the supercell in charge 
neutral. With the lattice constants fixed at the experimental 
values of a=4.138 A and c=28.64 A, the internal coordinates 
of all the atoms are fully relaxed until the residual forces on 
each atom are less than 0.02 eV/A. All these parameters en- 
sure the numerical errors in the calculated defect formation 
energies are less than 20 meV. 

Extensive research on the traditional diluted magnetic semi- 
conductors has indicated that their magnetic properties de- 
pend sensitively on the positions of magnetic impurities, 
namely, substitutional versus interstitial sites, in the host semi- 
conductors [17]. Therefore, we first examine the site prefer- 
ence of a serial of single 3d TM impurities (V, Cr, Mn, and 
Fe) in bulk Bi2Se3. To this end, we place the TM impuri- 
ties at all possible interstitial and substitutional sites in bulk 
Bi2Se3 and compare their formation energies, which can be 
computed using the following expression [38] 



AHf (TM) = E tot (TM) - E tot (bulk) 



(1) 



where E tot (TM) is the total energies of a supercell contain- 
ing one impurity, E to t(pulk) is the total energy of the equiv- 
alent supercell containing only host atoms, /i^ denotes the 
chemical potential for species i (host atoms or dopants), and 
rti indicates the corresponding number that have been added 
to or removed from the supercell. Given the host material 
being a binary compound, the values of ^ are subject to 
the following relations: (i) + 3/is- e = AHf(Bi2Ses), 
where AHf(Bi2Ses) is the formation energy of bulk Bi 2 Se3, 
to maintain equilibrium growth of Bi2Se3; (ii) jiBi ^ H°Bi> 
Use ^ A^e' anc ^ Mtm ^ Mtm to av °id precipitation of Bi, 
Se, and TM elementary substances, where n Bi , /i^ e , and fi^ M 
refer to the chemical potentials of the stable Bi, Se, and TMs 
crystals, respectively. We chose rhombohedral Bi, hexagonal 
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FIG. 2: (color online). Calculated formation energies of the most 
favorable charge states of single V, Cr, Mn, and Fe impurities doped 
at Bi substitutional sites in Bi2Se3 as a function of the Fermi level 
under (a) Bi-rich and (b) Se-rich conditions. 



Se, nonmagnetic body-centered-cubic (bcc) V, antiferromag- 
netic bcc Cr, antiferromagnetic face-centered-cubic Mn, and 
ferromagnetic bcc Fe as reference to evaluate the chemical po- 
tentials of the elements; (in) xfiTM+VI^Se ^ AHf(TM x Se y ) 
to ensure that all possible competing phases TM^Se^ can not 
precipitate. 

Our calculated formation energies indicate that Bi substitu- 
tional sites are strongly preferred by all the TM impurities we 
investigated, compared to Se substitutional sites and all possi- 
ble interstitial sites in bulk Bi 2 Se3. Interstitial sites in the van 
der Waals gap are the most favorable interstitial sites, while 
their formation energies are always larger than the dopants at 
substitutional Bi by at least 0.29 and 1.09 eV at Bi-rich and 
Se-rich conditions, respectively. The reported experimental 
observations of interstitial dopants in the van der Waals gap 
may be due to certain kinetic constraint, that is, the growth 
temperatures are not high enough for the dopants to overcome 
the energy barriers to incorporate into the substitutional sites. 
Fig. 1(b) shows the formation energies of the TM impuri- 
ties at Bi substitutional sites as a function of the chemical 
potentials of Bi and Se. The order of the formation ener- 
gies, AH f (V) < AH f (Cr) < AH f (Mn) < AH f (Fe), 
can be attributed to the size effect of the elements, that is, the 
larger the atomic radius of the impurity, the closer to that of 
the replaced Bi atom and thus the lower the formation energy. 
Among the four TM elements, V and Cr possess negative for- 
mation energies in the whole range of the accessible host ele- 
ment chemical potentials, suggesting that these two elements 
can be spontaneously introduced in bulk Bi2Se3. 

To further identify the optimal growth conditions, we inves- 
tigate the competition of the TM doping with the formation of 
intrinsic defects. The potential intrinsic defects include bis- 
muth vacancies Vbu selenium vacancies Vs e > an d antisite de- 
fects Bise and Sesi in bulk Bi2Se3. The most favorable in- 
trinsic defects and their formation energies are also shown in 
Fig. 1(b). We find that at extremely Bi-rich growth conditions, 
antisite defect Bis e i is favorable, and Bi2Se3 is expected to 
be intrinsic p type. However, at the other growth conditions, 
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FIG. 3: (color online), (a) Band structure of Bi2Se3 calculated for a 2 x 2 x 1 relaxed supercell with SOC. (b) Brillouin zone and high 
symmetry points of the 2 x 2 x 1 supercell. (c)-(j) Band structures of Cr doped and Fe doped Bi2Se3 calculated with: (c) and (d) structural 
relaxation but without SOC; (e) and (f) SOC but without structural relaxation; (g) and (h) SOC and structural relaxation; (i) and (j) SOC, 
structural relaxation, and GGA+U (U=3 eV). Structural relaxations were all performed without the inclusion of SOC. Projected partial density 
of states (DOS) on the dopants and the neighboring Se atoms are also shown in (e)-(h). 



selenium vacancies Vs e i an d antisite defects Scbu with both 
acting as donors, are dominant, rendering Bi 2 Se3 intrinsic n 
type. This result provides a good explanation to experimen- 
tal study [5, 6] and guidelines for growth control. Our results 
also indicate that the formation energies of V, Cr, and Mn are 
always lower than those of the intrinsic defects, whereas se- 
lenium vacancies (Ysei) an d antisite defects (Bis e i) are en- 
ergetically more favorable at Bi-rich conditions in the case 
of Fe doping. Therefore, in order to effectively dope Fe in 
bulk Bi 2 Se3, it is required to grow the materials at Se-rich 
conditions. However, the thermodynamic solubility of Fe in 
Bi 2 Se3 is expected to be limited, where the maximum doping 
concentration of Bi^^Fe^Ses is estimated to be x = 2.24% 
under Se-rich conditions at the melting point of Bi 2 Se3 of 979 
K, according to Ref. [38]. Experimentally, the effective dop- 
ing concentration of Fe in Bi 2 Se3 is reported to be less than 
2% [18], which is consistent with our result. 

In the above calculations, the whole systems are kept charge 



neutral, and the TM impurities, therefore, hold a 3+ charge 
state as they are introduced in Bi sites in Bi 2 Se3. Next, we 
study the formation energies of the doped TM impurities in 
various charge states. In this case, the definition of the impu- 
rity formation energy is modified as [38] 

AH f (TM, q) = E tot (TM, q) - E tot (bulk) - J^rc^ 

i 

+ ^ F +^ + AV), (2) 

where E tot (TM : q) is the total energy of a supercell contain- 
ing one impurity in charge state q and E F is the Fermi level 
with respect to the valence band maximum of the intrinsic 
bulk Bi 2 Se3 (Ey). AV is a potential alignment due to dif- 
ferent energy reference among the defect supercell and pure 
supercell. Figs. 2(a) and 2(b) show the calculated formation 
energies of the TM impurities in their favorable charge states 
as a function of Ep at Bi-rich and Se-rich growth conditions, 
respectively. The 3+ charge state of the TM impurities are 
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preserved when E F is smaller than 0.27, 0.28, 0.19, and 0.19 
eV for V, Cr, Mn, and Fe doped cases, respectively. At ex- 
treme n-type conditions, i.e., Ep approaching the conduction 
band minimal, the charge state of the TM impurities changes 
to 2+, which may result in the change of carrier concentra- 
tion in the materials when doping TM atoms. From Figs. 2(a) 
and 2(b), we obtain that the Se-rich condition is more optimal 
than the Bi-rich condition for doping TM in Bi2Se3. A sys- 
tematic investigation of the molecular beam epitaxy growth 
for Bi 2 Se3 films also indicates that the optimal growth can be 
achieved under a Se-rich atmosphere [39], which agrees with 
our results. 

We next focus our attention on the electronic structures of 
the TM doped Bi 2 Se3. Structural relaxations were all per- 
formed without the inclusion of SOC. The calculated band 
structures do not show any noticeable difference, in compari- 
son with the calculations with SOC included in the relaxation 
for Cr and Fe doped Bi2Se3, and the band gap differences are 
less than 0.01 eV. Figure 3(a) shows the calculated band struc- 
ture of pure Bi2Se3 with SOC, which gives a band gap of 0.32 
eV. In the TM doped cases, even without SOC, Cr, and Fe 
doping results in insulating magnetic states with energy band 
gaps of 0.28 [Fig. 3(c)] and 0.18 eV [Fig. 3(d)], respectively, 
whereas V and Mn create impurity states in the band gap of 
Bi 2 Se3, which completely close the band gap and make the 
materials metallic. With SOC, we obtain band gaps of 0.15 
and 0.03 eV for Cr and Fe doped Bi 2 Se3 as the dopants are 
fixed at the Bi sites without any structural relaxation [Figs. 
3(e) and 3(f)]. However, after the structural relaxation, the 
band gaps are further reduced, as shown in Figs. 3(g) and 3(h), 
to 0.010 and 0.028 eV for Cr and Fe doped Bi 2 Se 3 , respec- 
tively. This seemingly counterintuitive result can be attributed 
to structural relaxation induced enhancement of the electronic 
state hybridization between the TM dopants and their neigh- 
boring Se atoms. When all the atoms are fixed at their original 
crystal positions, due to the relatively small covalent radius of 
Cr and Fe compared with Bi, Cr and Fe states have little over- 
lap with the neighboring Se atoms, resulting in atomiclike lo- 
calized electronic states, as evidenced by the presence of the 
fairly flat impurity bands in the gap region shown in Figs. 3(e) 
and 3(f). After the structural relaxation, the distance between 
Cr (Fe) and the neighboring Se is decreased by about 0.4 A 
(0.3 A) and the impurity bands are substantially broadened, 
suggesting strong hybridization between the TM dopants and 
the neighboring Se atoms. A more detailed analysis on the 
projected partial density of states reveals that the d states of 
Cr and Fe hybridize the p states of Se. 

We further performed GGA+U calculations with U on the 
TM impurities ranging from 3 to 6 eV and J=0.87 eV to 
describe the strong electron-electron correlation in partially 
filled 3d TM elements. The results, as shown in Figs. 3(i) and 
3(j), indicate that the band gaps of Cr and Fe doped Bi 2 Se3 
increases to 0.025 and 0.028 eV for U=3 eV and 0.026 and 
0.028 eV for U=6 eV, respectively. Thus, electron-electron 
correlation could slightly enlarge the band gaps. 

Finally, we discuss the magnetic properties of the TM 



doped Bi 2 Se3. The magnetic ground states for single impu- 
rity doped Bi 2 Se3 have been calculated with the inclusion of 
SOC, the results indicate that both Cr and Fe prefer the di- 
rection perpendicular to the Bi 2 Se3 QLs, with the magnetic 
moments close to 3 and 5 respectively. Our results are 
different than the recent reports of in-plane magnetization 
easy axis in Ref. [40], where Fe adatoms on the surface of 
Bi 2 Se3(lll) and the intergrown composite crystal of FeySeg 
and Bi 2 Se3 [41], in which magnetic anisotropy is attributed 
to the intrinsic properties of ferromagnet FerSeg. The mag- 
netic coupling between the TM dopants is estimated by plac- 
ing two TM dopants in the supercell at various separation 
distances and calculating the total energy difference between 
antiferromagnetic (AFM) and ferromagnetic (FM) states of 
the dopants at the same distance. Our preliminary results re- 
veal that the magnetic coupling dominantly favors FM in the 
Cr doped Bi 2 Se3 while exclusively favors weak AFM in the 
Fe doped case. In the Cr doped case, the coupling strength 
[(Eafm — Efm)/ 2] is on the order of 10 meV for the two Cr 
atoms within the same QL at the first three nearest neighbor- 
ing distances and substantially reduced to the order of 1 meV 
or even less for the two Cr atoms in different QLs. We notice 
that studies of the competition between weak Localization and 
weak antilocalization from Lu et al [42] indicated that a fer- 
rromagnetic order can be formed in Cr doped Bi 2 Se3, which 
was later verified experimentally with a nominal Cr doping 
level of 23% [20]. More recently, ferromagnetism was re- 
ported to be observed in Cr doped Bi 2 Se3 [29], which agrees 
well with our result. 

In the Fe doped case, AFM is favorable while its coupling 
strength is very weak (on the order of 1 meV), even though the 
distance of two Fe atoms is close to 9.84 A. Experimentally, 
as we previously discussed, the doping concentration of Fe in 
Bi 2 Se3 is found to be lower than 2% [18], where the average 
distance between two uniformly distributed Fe atoms is about 
11.79 A. The corresponding coupling strength is predicted to 
be even less than 1 meV in this case. Hence, AFM in Fe doped 
Bi 2 Se3 may not be obvious to detect, especially in the low 
doping concentration regime, which agrees with experimental 
observations of nonmagnetism [24, 30]. In addition, heavily 
Fe doped Bi 2 Se3 often exhibit FM behavior [28], since Fe 
often tends to form clusters or become phase separated instead 
of being uniformly distributed during growth. 

In summary, we have studied the effects of magnetic dop- 
ing of a series of 3d transition metal elements in the TI Bi 2 Se3 
using first-principles calculations. Our calculated formation 
energies indicate that Bi substitutional sites are strongly pre- 
ferred by all the TM impurities. By examining the energetics 
of the dopants in the host material and the resulting electronic 
structures and magnetic properties, we have found that Cr and 
Fe doping preserves the insulating nature of the host TI in 
the bulk but the intrinsic band gap is substantially reduced 
due to the strong hybridization between the electronic states 
of the dopants and the neighboring Se atoms. The further in- 
vestigation on the magnetic coupling between the dopants has 
suggested that Cr doped Bi 2 Se3 is likely to be FM, while Fe 
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doped material tends to be weak AFM. Only Cr doped Bi 2 Se3 
with both ferromagnetic and insulating properties is promis- 
ing for realizing the quantized anomalous Hall effect. 

J.Z. and Y. Y. were supported by the MOST Project of China 
(Grant No. 2011CBA00100), NSF of China (Grants No. 
10974231, No. 11174337, and No. 11225418). W.Z. and 
D.X. were supported by the U.S. Department of Energy, Of- 
fice of Basic Energy Sciences, Materials Sciences and Engi- 
neering Division. W.Z. was also supported by the National 
Natural Science Foundation of China (grant No. 11034006). 
The calculations were performed at the Supercomputing Cen- 
ter of Chinese Academy of Sciences. 



* Electronic address: wzhu3@utk.edu 
t Electronic address: ygyao@bit.edu.cn 
[1] B. A. Bernevig and S.-C. Zhang, Phys. Rev. Lett. 96, 106802 

(2006) . 

[2] J. E. Moore and L. Balents, Phys. Rev. B 75, 121306 (2007). 
[3] L. Fu, C. L. Kane, and E. J. Mele, Phys. Rev. Lett. 98, 106803 

(2007) . 

[4] H. Zhang, C.-X. Liu, X.-L. Qi, X. Dai, Z. Fang, and S.-C. 

Zhang, Nat. Phys. 5, 438 (2009). 
[5] Y. Xia, D. Qian, D. Hsieh, L. Wray, A. Pal, H. Lin, A. Bansil, 

D. Grauer, Y. S. Hor, R. J. Cava, and M. Z. Hasan, Nat. Phys. 

5, 398 (2009). 

[6] D. Hsieh, Y Xia, D. Qian, L. Wray, J. H. Dil, F. Meier, J. Oster- 
walder, L. Patthey, J. G. Checkelsky, N. P. Ong, A. V. Fedorov, 
H. Lin, A. Bansil, D. Grauer, Y S. Hor, R. J. Cava, and M. Z. 
Hasan, Nature (London) 460, 1101 (2009). 

[7] M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. 82, 3045 (2010). 

[8] X.-L. Qi and S.-C. Zhang, Rev. Mod. Phys. 83, 1057 (2011). 

[9] Q. Liu, C.-X. Liu, C. Xu, X.-L. Qi, and S.-C. Zhang, Phys. Rev. 
Lett. 102, 156603 (2009). 
[10] X.-L. Qi, R. Li, J. Zang, and S.-C. Zhang, Science 323, 61 
(2009). 

[11] I. Garate, and M. Franz, Phys. Rev. Lett. 104, 146802 (2010). 
[12] L. A. Wray, S.-Y Xu, Y Xia, D. Hsieh, A. V. Fedorov, Y S. 

Hor, R. J. Cava, A. Bansil, H. Lin, and M. Z. Hasan, Nat. Phys. 

7,32(2011). 

[13] Y L. Chen, J.-H. Chu, J. G. Analytis, Z. K. Liu, K. Igarashi, 
H.-H. Kuo, X. L. Qi, S. K. Mo, R. G. Moore, D. H. Lu, M. 
Hashimoto, T. Sasagawa, S. C Zhang, I. R. Fisher, Z. Hussain, 
and Z. X. Shen, Science 329, 659 (2010). 

[14] R. Yu, W. Zhang, H.-J. Zhang, S.-C Zhang, X. Dai, and Z. 
Fang, Science 329, 61 (2010). 

[15] C.-Z. Chang, J.-S. Zhang, M.-H. Liu, Z.-C. Zhang, X. Feng, 
K. Li, L.-L. Wang, X. Chen, X. Dai, Z. Fang, X.-L. Qi, 
S.-C. Zhang, Y Wang, K. He, X.-C. Ma, and Q.-K. Xue, 
arXiv:1108.4754vl. 

[16] Z. Qiao, S. A. Yang, W Feng, W-K. Tse, J. Ding, Y Yao, J. 
Wang, and Q. Niu, Phys. Rev. B 82, 161414(R) (2010); J. Ding, 
Z. Qiao, W. Feng, Y Yao, and Q. Niu, ibid. 84, 195444 (2011). 

[17] T. Jungwirth, J. Sinova, J. Masek, J. Kucera, and A. H. Mac- 
Donald, Rev. Mod. Phys. 78, 809 (2006). 

[18] J. J. Cha, J. R. Williams, D. Kong, S. Meister, H. Peng, A. 
J. Bestwick, P. Gallagher, D. Goldhaber-Gordon, and Y Cui, 
Nano Lett. 10, 1076 (2010). 



[19] H. Jin, J. Im, and A. J. Freeman, Phys. Rev. B 84, 134408 
(2011). 

[20] M. Liu, J. Zhang, C.-Z. Chang, Z. Zhang, X. Feng, K. Li, K. 
He, L.-l. Wang, X. Chen, X. Dai, Z. Fang, Q.-K. Xue, X. Ma, 
and Y. Wang, Phys. Rev. Lett. 108, 036805 (2012). 

[21] J. S. Dyck, P. Hajek, P. Lost'ak, and C. Uher, Phys. Rev. B 65, 
115212 (2002); J. S. Dyck, C. Drasar, P. Lost'ak, and C. Uher, 
ibid. 71, 115214 (2005); Z. Zhou, Y-J. Chien, and C. Uher, 
ibid. 74, 224418 (2006); Y-J. Chien, Ph.D. thesis, University 
of Michigan, 2007. 

[22] J. Choi, H.-W Lee, B.-S. Kim, S. Choi, J. Choi, J. H. Song, and 
S. Cho, J. Appl. Phys. 97, 10D324 (2005). 

[23] P. Larson and W R. L. Lambrecht, Phys. Rev. B 78, 195207 
(2008). 

[24] V. A. Kulbachinskii, A. Y Kaminskii, K. Kindo, Y Narumi, 
K. Suga, P. Lostak, and P. Svanda, Physica (Amsterdam) 31 IB, 
292 (2002). 

[25] J. Choi, S. Choi, J. Choi, Y Park, H.-M. Park, H.-W Lee, B.-C. 

Woo, and S. Cho, Phys. Status Solidi B 241,1541 (2004). 
[26] Y S. Hor, P. Roushan, H. Beidenkopf, J. Seo, D. Qu, J. G. 

Checkelsky, L. A. Wray, D. Hsieh, Y Xia, S.-Y Xu, D. Qian, 

M. Z. Hasan, N. P. Ong, A. Yazdani, and R. J. Cava, Phys. Rev. 

B 81, 195203 (2010). 
[27] C. Niu, Y Dai, M. Guo, W Wei, Y Ma, and B. Huang, Appl. 

Phys. Lett. 98, 252502 (2011). 
[28] Y H. Choi, N. H. Jo, K. J. Lee, J. B. Yoon, C. Y You, and M. 

H. Jung, J. Appl. Phys. 109, 07E312 (2011). 
[29] P. P. J. Haazen, J.-B. Laloe, T. J. Nummy, H. J. M. Swagten, 

P. Jarillo-Herrero, D. Heiman, and J. S. Moodera, Appl. Phys. 

Lett. 100, 082404 (2012). 
[30] Y. Sugama, T. Hayashi, H. Nakagawa, N. Miura, and V. A. Kul- 

bachnskii, Physica (Amsterdam) 298B, 531 (2001). 
[31] Z. Salman, E. Pomjakushina, V. Pomjakushin, A. Kanigel, K. 

Chashka, K. Conder, E. Morenzoni, T. Prokscha, K. Sedlak, and 

A. Suter, arXiv:1203.4850vl. 

[32] P. E. Blochl, Phys. Rev. B 50, 17953 (1994). 
[33] J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 
3865 (1996). 

[34] G. Kresse and J. Hafner, Phys. Rev. B 48, 13115 (1993). 

[35] G. Kresse and J. Furthmiiller, Phys. Rev. B 54, 11169 (1996). 

[36] Due to the important relativistic effects on Bi2Se3, SOC was 
taken into account in our calculations. The test value of forma- 
tion energy for Bi 2 Se3 with SOC is about 17% different from 
that without SOC. 

[37] S. B. Zhang and J. E. Northrup, Phys. Rev. Lett. 67, 2339 
(1991). 

[38] C. G. Van de Walle and J. Neugebauer, J. Appl. Phys. 95, 3851 
(2004). 

[39] C.-L. Song, Y-L. Wang, Y-P. Jiang, Y Zhang, C.-Z. Chang, 
L. Wang, K. He, X. Chen, J.-F. Jia, Y Wang, Z. Fang, X. Dai, 
X.-C. Xie, X.-L. Qi, S.-C. Zhang, Q.-K. Xue, andX. Ma, Appl. 
Phys. Lett. 97, 143118 (2010). 

[40] J. Honolka, A. A. Khajetoorians, V. Sessi, T. O. Wehling, S. 
Stepanow, J.-L. Mi, B. B. Iversen, T. Schlenk, J. Wiebe, N. 

B. Brookes, A. I. Lichtenstein, Ph. Hofmann, K. Kern, and R. 
Wiesendanger, Phys. Rev. Lett. 108, 256811 (2012). 

[41] H. Ji, J. M. Allred, N. Ni, J. Tao, M. Neupane, A. Wray, S. Xu, 
M. Z. Hasan, and R. J. Cava, Phys. Rev. B 85, 165313 (2012). 

[42] H.-Z. Lu, J. Shi, and S.-Q. Shen, Phys. Rev. Lett. 107, 076801 
(2011). 



